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Abstract 

We report the results of the femtoscopic analysis of pairs of identical pions measured in p-Pb col¬ 
lisions at i/snn = 5.02 TeV. Femtoscopic radii are determined as a function of event multiplicity 
and pair momentum in three spatial dimensions. As in the pp collision system, the analysis is com¬ 
plicated by the presence of sizable background correlation structures in addition to the femtoscopic 
signal. The radii increase with event multiplicity and decrease with pair transverse momentum. When 
taken at comparable multiplicity, the radii measured in p-Pb collisions, at high multiplicity and low 
pair transverse momentum, are 10-20% higher than those observed in pp collisions but below those 
observed in A-A collisions. The results are compared to hydrodynamic predictions at large event 
multiplicity as well as discussed in the context of calculations based on gluon saturation. 
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1 Introduction 

The Large Hadron Collider (LHC) OJ] delivered Pb-Pb collisions at ^/snn = 2.76 TeV in 2010 and in 
2011. Several signatures of a quark-gluon plasma were observed, including a strong suppression of high- 
p-\ particle production (”jet quenching”) [2 4|. as well as collective behavior at low pj [ 5 , i§] which is 
well described by hydrodynamic models with a low shear viscosity to entropy ratio. A comparison to 
reference results from pp collisions at \fs = 0.9, 2.76, and 7 TeV shows that these effects cannot be 
described by an incoherent superposition of nucleon-nucleon collisions. As such, they can be interpreted 
as final-state phenomena, characteristic of the new state of matter lEffl] created in heavy-ion collisions. 
In order to verify the creation of such a state, p-Pb collisions at ^/.vnn = 5.02 TeV, where a quark-gluon 
plasma is not expected to form, were provided by the LHC. In particular, cold nuclear matter effects, such 
as gluon saturation, which are expected to influence a number of observables, are being investigated [11]. 

One of the observables characterizing the bulk collective system is the size of the particle-emitting region 
at freeze-out which can be extracted with femtoscopic techniques IU 2j,[i3|]. In particular, two-particle cor¬ 
relations of identical pions (referred to as Bose-Einstein, or Hanbury-Brown Twiss ”HBT”, correlations) 
provide a detailed picture of the system size and its dependence on the pair transverse momentum and the 
event multiplicity. Femtoscopy measures the apparent width of the distribution of relative separation of 
emission points, which is conventionally called the “radius parameter.” The radius can be determined in¬ 
dependently for three directions: long along the beam axis, out along the pair transverse momentum, and 
side, perpendicular to the other two. Such measurements were performed at the LHC for central Pb-Pb 
collisions [ 14] as well as for pp collisions at yfs = 0.9 and 7 TeV [15-18], and compared to results from 
heavy-ion collisions at lower collision energies. Two clear trends were found: (i) In A-A collisions all 
radii scale approximately linearly with the cube root of the final state charged particle multiplicity den¬ 
sity at midrapidity (dV c h/dq) 1/3 for all three radii separately, consistently with previous findings |JL3]. 
For pp collisions, the radii scale linearly with the cube root of charged particle multiplicity density as 
well, however the slope and intercept of the scaling line are clearly different than for A-A. (ii) A signif¬ 
icant, universal decrease of the radii with pair momentum has been observed in A-A collisions, while 
the analogous trend in pp depends on the considered direction ( out , side, or long ) and event multiplicity. 
A transverse momentum dependence of the radii similar to A-A was observed for the asymmetric d-Au 
collision system at RHIC [ 19i, & The one-dimensional radii extracted from the ALICE 3-pion cumu- 
lant analysis were also investigated in pp, p-Pb, and Pb-Pb collisions at the LHC [21]. For the p-Pb 
system, at a given multiplicity, the radii were found to be 5-15% larger than those in pp, while the radii 
in Pb-Pb were 35-55% larger than those in p-Pb. 

The A-A pion femtoscopy results are interpreted within the hydrodynamic framework as a signature 
of collective radial flow. Models including this effect are able to reproduce the ALICE data for central 
collisions [22j,|23]. Attempts to describe the pp data in the same framework have not been successful so 
far and it is speculated that additional effects related to the uncertainty principle may play a role in such 
small systems [24]. In p-A collisions, hydrodynamic models which assume the creation of a hot and 
dense system expanding hydrodynamically predict system sizes larger than those observed in pp, and 
comparable to those observed in lower-energy A-A collisions at the same multiplicity [24, 25]- However 
such models have an inherent uncertainty of the initial state shape and size, which can also differ between 
pp and peripheral A-A collisions. 

Alternatively, a model based on gluon saturation suggests that the initial system size in p-A collisions 
should be similar to that observed in pp collisions, at least in the transverse direction \2Q, 27]. At that 
stage both systems are treated in the same manner in the Color Glass Condensate model (CGC), so that 
their subsequent evolution should lead to comparable radius parameter at freeze-out. Ref. [281] suggests a 
(small) Yang-Mills evolution in addition. The observation of a larger size in the p-A system with respect 
to pp would mean that a comparable initial state evolves differently in the two cases, which is not easily 
explained within the CGC approach alone. The d-Au data measured at RHIC suggest that hydrodynamic 
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evolution may be present in such a system, while the ALICE 3-pion analysis at the LHC [21] leaves room 
for different interpretations. The pion femtoscopic radii as a function of pair transverse momentum from 
p-Pb collisions at v Ann = 5.02 TeV, which are reported in this paper, provide additional constraints on 
the validity of both approaches. 


The paper is organized as follows: In Sec. [2] the data-taking conditions, together with event and track 
selections are described. The femtoscopic correlation function analysis, as well as the extraction of the 
radii and associated systematic uncertainties, and the discussion of the fitting procedure are explained 
in Sec. [3] In Sec. [J] the results for the radii are shown and compared to model predictions. Section [5] 
concludes the paper. 


2 Data taking and track reconstruction 


The LHC delivered p-Pb collisions at the beginning of 2013 at ^Ann = 5.02 TeV (4 TeV and 1.58 TeV 
per nucleon for the p and Pb beams, respectively). The nucleon-nucleon center-of-mass system is shifted 
with respect to the ALICE laboratory system by 0.465 unit of rapidity in the direction of the proton beam. 


The ALICE detector and its performance are described in Refs. 1291.1301. The main triggering detector 
is the V0, consisting of two arrays of 32 scintillator counters, which are installed on each side of the 
interaction point and cover 2.8 < r|i a b <5.1 (V0A, located on the Pb-remnant side), and —3.7 < rp a b < 
— 1.7 (V0C). The minimum-bias trigger requires a signal in both V0 detectors within a time window 
that is consistent with the collision occurring at the center of the ALICE detector. Additionally, specific 
selection criteria to remove pile-up collisions are applied [300. Approximately 80 million minimum-bias 
events were analyzed. 


The analysis was performed in multiplicity classes, which were determined based on the signal from 
the V0A detector, located along the Pb-going beam. This ensures that the multiplicity determination 
procedure uses particles at rapidities significantly different from the ones used for the pion correlation 
analysis, avoiding potential auto-correlation effects. Events are grouped in four multiplicity classes: 
0-20%, 20-40%, 40-60%, and 60-90%, defined as fractions of the analyzed event sample sorted by 
decreasing V0A signal, which is proportional to the multiplicity within the acceptance of this detector. 
Table [Qshows the multiplicity class definitions and the corresponding mean charged-particle multiplicity 
densities (d/V c h/dr|) averaged over |r|i a b| < 0.5 as obtained using the method presented in Ref. 113111 . The 
(dV c h/dri) values are not corrected for trigger and vertex-reconstruction inefficiency, which is about 4% 
for non-single diffractive events [31]. 


Event class 

(dV C h/dri) 

Ifiiabl < 0.5,p T > 0 GeV/c 

60-90% 

8.2±0.3 

40-60% 

16.1 ±0.4 

20-40%: 

23.2 ±0.5 

0-20% 

35.5 ±0.8 


Table 1: Definition of the V0A multiplicity classes as fractions of the analyzed event sample and their correspond¬ 
ing (dV c h/dr|(|r|i a b| < 0.5.p-\ > 0)). The given uncertainties are systematic only since the statistical uncertainties 
are negligible. 

Charged track reconstruction is performed using the Time Projection Chamber (TPC) and the Inner 
Tracking System (ITS). The TPC is a large volume cylindrical gaseous tracking chamber, providing 
information of particle trajectories and their specific energy loss. The read-out chambers mounted on 
the endcaps are arranged in 18 sectors on each side (covering the full azimuthal angle) measuring up to 
159 samples per track. The TPC covers an acceptance of |r|i a b| <0.8 for tracks which reach the outer 


3 












Two-pion femtoscopy in p-Pb collisions 


ALICE Collaboration 


radius of the TPC and |r|i a b| <1.5 for shorter tracks. The ITS is composed of position-sensitive silicon 
detectors. It consists of six cylindrical layers: two layers of Silicon Pixel Detector (SPD) closest to the 
beam pipe covering |r|i a b| < 2.0 and |r|i a b| <1.4 for inner and outer layers respectively, two layers of 
Silicon Drift Detector in the middle covering |r|i a b| < 0.9, and two layers of Silicon Strip Detector on 
the outside covering |r|i a b| < 1.0. The information from the ITS is used for tracking and primary particle 
selection. The momentum of each track is determined from its curvature in the uniform magnetic field 
of 0.5 T oriented along the beam axis, provided by the ALICE solenoidal magnet. 


The primary-vertex position is determined with tracks reconstructed in the ITS and TPC as described in 
Ref. 13211 . Events are selected if the vertex position along the beam direction is within ±10 cm of the 
center of the detector. This ensures a uniform acceptance in rp a b. 


Each track is required to exploit signals in both TPC and ITS. The track segments from both detectors 
have to match. Additionally, each track is required to have at least one hit in the SPD. A TPC track 
segment is reconstructed from space points (clusters). Each track is required to be composed of at least 
50 out of the 159 such clusters. The parameters of the track are determined by a Kalman fit to the set 
of TPC+ITS clusters. The quality of the fit % 2 was required to be better than 4 per cluster in the TPC 
and better than 36 in ITS. Tracks that show a kink topology in the TPC are rejected. To ensure that 
dominantly primary-particle tracks are selected, the distance of closest approach to the primary vertex 
is required to be closer than 2.0 cm in the longitudinal direction and (0.0105 + 0.0350 • /; T 11 ) cm, with 
pj in GeV/c, in the transverse direction. The kinematic range of particles selected for this analysis is 
0.12 < px < 4.0 GeV/c and |rji a b| < 0.8. 

The Time-Of-Flight (TOF) detector is used together with the TPC for pion identification. The TOF is a 
cylindrical detector of modular structure, consisting of 18 azimuthal sectors divided into 5 modules along 
the beam axis at a radius r ~ 380 cm. The active elements are Multigap Resistive Chambers (MRPC). 
For both TPC and TOF, the signal (specific energy loss AE/Ax for the TPC and the time-of-flight for 
TOF) for each reconstructed particle is compared with the one expected for a pion. The difference is 
confronted with the detector resolution. The allowed deviations vary between 2 to 5a for the TPC and 2 
to 3a for TOF depending on the momentum of the particle. The selection criteria are optimized to obtain 
a high-purity sample while maximizing efficiency, especially in the regions where the expected signal 
for other particles (electrons, kaons, and protons for the TPC, kaons for TOF) approaches the pion value. 
The purity of the pion sample is above 98%. 


The accepted particles from each event are combined to pairs. The two-particle detector acceptance 
effects, track splitting and track merging, are present. Track splitting occurs when a single trajectory 
is mistakenly reconstructed as two tracks. AFICE tracking algorithm has been specifically designed to 
suppress such cases. In a rare event when splitting happens, two tracks are reconstructed mostly from the 
same clusters in the AFICE TPC. Therefore pairs which share more than 5% of clusters are removed from 
the sample. Together with the anti-merging cut described below this eliminates the influence of the split 
pairs. Track merging can be understood as two-particle correlated efficiency and separation power. In 
the ALICE TPC, two tracks cannot be distinguished if their trajectories stay close to each other through 
a significant part of the TPC volume. Although this happens rarely such pairs by definition have low 
relative momentum and therefore their absence distorts the correlation function in the signal region. 
Track splitting and track merging are taken into account and corrected with the procedure described in 
Ref. [16]. The effect of the two-particle detector acceptance on the final results is similar to what was 
observed in pp and is limited to low pair relative momentum, where it slightly affects the shape of the 
correlation function. However, in p-Pb collisions the femtoscopic effect is an order of magnitude wider 
than any region affected by this inefficiency and, as a consequence, the extracted radii are not affected 
by the two-hack acceptance. 
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Fig. 1: (Color online) Projections of the three-dimensional n + n + correlation functions for three selected multi¬ 
plicity and h r ranges along the out (top), side (middle), and long (bottom) direction. The other components are 
integrated over the four bins closest to zero in their respective q directions. 


3 Correlation function analysis 

3.1 Construction of the correlation function 


The correlation function C(pi .pi) of two particles with momenta pi and p 2 is defined as 


C(P1,P2) 


^(PmPs) 

fl(Pl.P2)' 


( 1 ) 


The signal distribution A is constructed from pairs of particles from the same event. The background 
distribution B should be constructed from uncorrelated particles measured with the same single-particle 
acceptance. It is built using the event mixing method with the two particles coming from two different 
events for which the vertex positions in beam direction agree within 2 cm and the multiplicities differ by 
no more than 1 /4 of the width of the given event class. The denominator is normalized to the number of 
entries in the numerator, so that the absence of correlation gives a correlation function at unity. 


The femtoscopic correlation is measured as a function of the momentum difference of the pair q = P 2 — Pi 
as 


C(q) 


A(q) 

5(q)’ 


( 2 ) 
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q (GeV/c) 

Fig. 2: (Color online) First three non-vanishing components of the SH representation of the n + n + correlation 
functions for three multiplicity and kj ranges, l = 0, m = 0 (top), / = 2, m = 0 (middle), and 1 = 2, m = 2 (bottom). 


where the dependence on the pair total transverse momentum kj = |pi.t + P2,t|/2 is investigated by per¬ 
forming the analysis in the following ranges in kj: 0.2-0.3, 0.3-0.4, 0.4-0.5, 0.5-0.6, 0.6-0.7, 0.7-0.8, 
and 0.8-1.0 GeV/c. The kj ranges are the same for each multiplicity class resulting in 28 indepen¬ 
dent correlation functions overall. Systematic uncertainties on the correlation functions are discussed in 
Sec.[3gl 


The momentum difference q is evaluated in the Longitudinally Co-Moving System (LCMS) frame in 
which the total longitudinal pair momentum vanishes: pi l + P2.L = 0 . similarly to previous measure¬ 
ments in small systems lllal . In Fig. Q] correlation functions are shown, projected over 128 MeV/c-wide 
slices along the q oul , ^ s ide> and c/| ong axes. An enhancement at low relative momentum is seen in all pro¬ 
jections. The width of this correlation peak grows with decreasing multiplicity or with increasing kj. The 
femtoscopic effect is expected to disappear at large q = |q|, with the correlation function approaching 
unity. We observe, especially for large kj and small multiplicities, that the correlation function is not flat 
in this region and has different values in different projection^- The cause may be non-femtoscopic corre¬ 
lations, which are presumably also affecting the shape of the correlation function in the femtoscopic (low 
q) region. This issue is a major source of systematic uncertainty on the extracted radii and is discussed 


1 We note that the overall normalization of the correlation function is a single value for the full three-dimensional object and 
cannot be independently tuned in all projections. 
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q (GeV/c) 


Fig. 3: (Color online) Dependence of the SH components of the correlation function on event multiplicity in a 
broad relative momentum range. 


in detail in Sections [3721 and [331 


The pair distributions and the correlation function can be represented in spherical harmonics (SH) 1331 . 


34] alternatively to the traditionally-used Cartesian coordinates. All odd-/ and odd-m components of 


such a representation vanish for symmetry reasons. The important features of the correlation function 
are fully captured by the following ones: / = 0, m = 0 is sensitive to the overall size of the pion source, 
/ = 2, m = 0 is sensitive to the difference between the longitudinal and transverse sizes, and / = 2, m = 2 
reflects the difference between the sidewards and outwards transverse radii. Therefore, three independent 
sizes of the source can also be extracted from these three SH components. 


In Fig. [2] we show the first three non-vanishing components of the spherical harmonics representation 
corresponding to the correlation functions shown in Fig. [H In the (0,0) component the enhancement 
at low-g is clearly visible, decreasing (increasing) in width with multiplicity (kj). The other two com¬ 
ponents, (2,0) and (2,2), also show structures in this region, indicating that the source shape is not 
spherically symmetric in the LCMS frame. 


3.2 Non-femtoscopic structures 

As mentioned in the discussion of Figs. Q] and |2j a significant non-femtoscopic correlation is observed 
in the range in q that is much larger than the characteristic width of the femtoscopic effect. As an ex¬ 
ample, in Fig. [3] we show the correlation in the SH representation up to 2.0 GeV/c in q. For the lowest 
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Fig. 4: (Color online) First three non-vanishing components of the SH representation of the 7t + 7t + correlation 
functions for a selected multiplicity and kj range, compared to a calculation from EPOS 3.076 1351, 1X61 (generator 
level only) and PYTHIA 6.4 Perugia-0 tune |37| '38 ] lor pp at i Js = 7 TeV (full simulation with detector response). 


multiplicity, and to a smaller degree at higher multiplicities, a significant slope in low q region is seen 
in the (0,0) component and a deviation from zero in the (2,0) component up to approximately 1 GeV/c. 
Similar correlations have been observed by ALICE in pp collisions & They were interpreted, based 
on Monte-Carlo model simulations, to be a manifestation of mini-jets, the collimated fragmentation of 
partons scattered with modest momentum transfer. The lowest multiplicities observed in p-Pb colli¬ 
sions are comparable to those in pp collisions at y/s = 7 TeV. Therefore a similar interpretation of the 
non-femtoscopic correlations in this analysis is natural. Similar structures have been observed in d-Au 
collisions by STAR lfl9h . This picture is corroborated by the analysis using the 3-pion cumulants, where 
expectedly the mini-jet contribution is suppressed [21]. 


Two important features of the non-femtoscopic correlation affect the interpretation of our results. First, 
it is a broad structure, extending up to 1 GeV/c and we have to assume that it also extends to 0 GeV/c 
in q. Therefore it affects the extracted femtoscopic radii and has to be taken into account in the fitting 
procedure. It can be quantified in the large q region and then extrapolated, with some assumptions, to the 
low q region, under the femtoscopic peak. The procedure leads to a systematic uncertainty. Secondly, 
it becomes visibly larger as multiplicity decreases and also as kj increases, which is consistent with the 
mini-jet-like correlation. 


The background was studied in Monte Carlo models, such as, for p-Pb collisions, AMPT [39], HIJING 
[ 0 ], DPMJET 0 , EPOS 3.076 00 ], and PYTHIA 6.4 (Perugia-0 tune) [ 3,0 for pp collisions 
at similar multiplicities. In all cases the Monte Carlo correlation functions exhibit significant structures 
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similar to the long-range effects observed in data, which is another argument for their non-femtoscopic 
origin. However, quantitative differences in the magnitude and shape of these structures when compared 
to those observed in data are seen for AMPT, HIJING and DPMJET. These models are therefore un¬ 
suitable for a precise characterization of the background, which is needed for the fitting procedure. The 
only models that qualitatively describe the features of the background (enhancement at low q, growing 
with kj, and falling with multiplicity) are EPOS 3.076 and PYTHIA 6.4 (Perugia-0 tune), which was also 
used in the pp analysis G3]. We note that PYTHIA simulation included full detector response modelling, 
while the EPOS 3.076 one did not. The comparison with data is shown in Fig. [4] The behavior of the 
correlation is well reproduced above 0.5 GeV/c in q, where non-femtoscopic correlations are expected 
to dominate. At low q, below 0.3 GeV/c, the data and models diverge, which is expected, as the femto- 
scopic correlations are not included in the model calculation. Above 0.3 GeV/c, EPOS reproduces the 
(0,0) component well, PYTHIA slightly overestimates the data. For the (2,0) and (2,2) components, 
which describe the three-dimensional shape of the non-femtoscopic correlations, PYTHIA is closer to 
the experimental data. Overall, for like-sign pairs, both models are reasonable approximations of the 
non-femtoscopic background. We use these models to fix the background parameters in the fitting pro¬ 
cedure. 

Similarly to the pp analysis Q, the unlike-sign pairs have also been studied. We found that correlations 
in the (0,0) component of PYTHIA are slightly larger than in data in the femtoscopic region for all kj 
ranges, and similar to data in the (2,0) and (2,2) components. EPOS was found to reasonably describe 
the unlike-sign pairs for low kj ranges and has smaller correlation than data in the (0,0) component in 
higher kj ranges. 


3.3 Fitting the correlation functions 

The space-time characteristics of the source are reflected in the correlation function 

C(q) = J S(r,q)\W(q,r)\ 2 d\ (3) 


where r is the pair space-time separation four-vector. S is the source emission function, interpreted as 
a probability density describing the emission of a pair of particles with a given relative momentum and 
space-time separation. V P is the two-particle interaction kernel. 


Previous femtoscopy studies in heavy-ion collisions at SPS [42], RHIC 143-4^1 and at the LHC 114] 
used a Gaussian static source S 


S{r) 


exp 


r 2 r 2 

out _ side 

4/? g 2 4R g 2 
n/r out ^'Lide 



(4) 


The R G ut , P G dc , and Rf ona parameters describe the single-particle source size in the LCMS in the out, side, 
and long directions, respectively. 


The Gaussian source provides a commonly used approximation of the source size and was used to com¬ 
pare to other experimental results, especially the ones coming from A-A collisions, where the source 
shape is more Gaussian than in small systems. While pursuing the standard procedure with the Gaussian 
assumption, we also carefully look for any deviations between the fit function and data which might 
suggest a significantly non-Gaussian shape of the source, which would be an important similarity to the 
pp case. 


In the analysis of pp collisions by ALICE a Gaussian is used together with other source shapes, 
exponential and Lorentzian 111611 . A Lorentzian parametrization in the out and long directions and a Gaus¬ 
sian parametrization in the side direction were found to fit the data best according to % 2 /ndf. Therefore, 
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we use this source parametrization also in the analysis of p-Pb collisions 


S(r) 


1 


'out 


+ R 


iE 2 

out 


exp 


' side 


1 


4R g 

4A side 


,+Rf. 


2 ' 


(5) 


long ' long 

The corresponding source sizes in out and long are R l mil and Rf ona , while for the side direction the size 
parameter R[’ lik is the same as in the Gaussian case. 


For identical pions, which are bosons, V P must be symmetrized. Since charged pions also interact via 
the Coulomb and strong Final State Interactions (FSI), I'Ll 2 corresponds to the Bethe-Salpeter ampli¬ 
tude & For like-sign pion pairs the contribution of the strong interaction is small for the expected 
source sizes 150], and is neglected here. The used l P therefore is a symmetrized Coulomb wave function. 
It is approximated by separating the Coulomb part and integrating it separately, following the procedure 
known as Bowler-Sinyukov fitting 1151 52], which was used previously for larger sizes observed in Pb- 
Pb [14] as well as smaller sizes observed in pp collisions 1611 . In this approximation the integration of 
Eq. © with S given by Eq. © results in the following functional form for the correlation function which 
is used to fit the data 


Cf(q) 


= (I-*) 

+ KKc(q) 



( 6 ) 


The function Kc(q ) is the Coulomb part of the two-pion wave function integrated over the spherical 
Gaussian source with a fixed radius. The value of this radius is chosen to be 2 fm. Its uncertainty has 
systematic effects on the final results (see Sec. 13.41) . This form of the correlation function from Eq. ([ 6 ]) is 
denoted in the following as GGG. Similarly for the source shape given by Eq. ([5]) the correlation function 
is 


C f (q) = (1-X) 


+ 


IK C 


1 +exp 




E 2 a 2 
long “ long 


(V) 


It has an exponential shape in out and long and a Gaussian shape in the side direction. Therefore it 
is referred to as EGE form of the correlation function. Parameter A in Eqs. © and © represents the 
correlation strength. 

Additionally, a component G describing non-femtoscopic correlations needs to be introduced. There is 
no a priori functional form which can be used for this component. Several of its features can be deduced 
from the correlations shown in Figs. [T} El and [3] it has to allow for different shapes in the out, side, and 
long directions; in the ( 0 , 0 ) component it has to extrapolate smoothly to low q and have a vanishing 
slope at q = 0. Since this structure is not known, maximum information about its shape and magnitude 
should be gained from an observation of the raw correlation functions and the corresponding effects in 
Monte Carlo, in as many representations as possible. It is therefore crucial to simultaneously use the 
Cartesian and SH representations as they provide complementary ways to study the correlation shape. 


An ad-hoc, Monte-Carlo driven parametrization of the non-femto background that reasonably describes 
the correlation function is G, composed of three independent 1 -dimensional functions G[J (Gaussian 
plus fixed constant), G 1 ] (Gaussian plus variable constant) and G 2 (Gaussian plus an additional linear 
component) 


G 8 ( 9 ) 
G ° 2 (q) 


Nq 


1 +a 0 exp - 


2K) 2 


0 $ • exp 


2 ( a «) 2 


+ P5, 


(8) 

(9) 
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Fig. 5: (Color online) First three non-vanishing components of the SH representation of the n + n + correlation 
functions for three multiplicity and kj combinations, l = 0 , m = 0 (top), I = 2 , m = 0 (middle), and l = 2 , m = 2 
(bottom). The lines show the corresponding components of the Gaussian (GGG) fit. 

^2(9) = a 2 ‘ ex P 2 (a 2 ) 2 ) “*"P2+'Y2‘‘7; (10) 

where Nq, a[j, a{], a 2 , o?, a 2 , <T 2 , and y 2 are fixed to the values obtained from fits to Monte Carlo events, 
separately for each multiplicity and kj range. In the fit procedure the |3 2 and P? parameters are kept free. 
This results in the following fit formula 

C(q)=A-C f (q)-[Q°(^)-To 0 (e,c P )+^(^)-T 2 0 (e,(p)+n2(^).T 2 2 (e,cp)], (11) 

where N is the overall normalization factor and T o °(0,(p), 7 2 °(0,cp), and L, 2 (0. (p) are the real parts of the 
relevant spherical harmonic functions. 

The fit is performed with the log-likelihood method in three dimensions for the Cartesian representation. 
The Gaussian fit reproduces the overall width of the femtoscopic correlation in all cases. The background 
component describes the behavior of the correlation at large q, but can also have non-zero correlation at 
0 in q. 

A corresponding fit is also performed for the SH representation of the correlation, which is shown in 
Fig. [5] The formula from Eq. © or Eq. (O (for the GGG fit or the EGE fit, respectively) is numerically 
integrated on a cp — 0 sphere for each q bin, with proper YJ n weights, to produce the three components 
of the SH decomposition. Statistical uncertainties on each component as well as the covariance matrix 
between them are taken into account in this simultaneous fit to the three histograms. The results are 
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Fig. 6: (Color online) First three non-vanishing components of the SH representation of the n + n + correlation 
functions for three multiplicity and kj combinations, l = 0, m = 0 (top), / = 2, m = 0 (middle), and l = 2, m = 2 
(bottom). The lines show the corresponding components of the EGE fit. 


shown in Fig. [5] The fit describes the general direction-averaged width of the correlation function, 
shown in the upper panel. The background component Q. describes the behavior at large q but also 
contributes to the correlation at low q. The shape in three-dimensional space, captured by the (2,0) and 
(2,2) components, is also a combination of the femtoscopic and non-femtoscopic correlations. 

Overall the GGG fit describes the width of the correlation but the data at low q are not perfectly re¬ 
produced, which can be attributed to the limitations of the Bowler-Sinyukov formula as well as to the 
non-Gaussian, long-range tails which are possibly present in the source. Some deviations from the pure 
Gaussian shape can also be seen for the long direction for the higher multiplicities. The EGE fit (Eq. ([7])) 
better reproduces the correlation peak in the (0,0) component, as shown in Fig. [6] The (2,0) and (2,2) 
components show similar quality of the fit. The % 2 values for both fits are comparable. 

3.4 Systematic uncertainties on the radii 

The analysis was performed separately for positively and negatively charged pions. For the practically 
zero-net-baryon-density system produced at the LHC they are expected to give consistent results. Both 
datasets are statistically consistent at the correlation function level. 

The main contributions to the systematic uncertainty are given in Table [2] for GGG radii and in Table [3] 
for EGE radii. 

We used two alternative representations (Cartesian and spherical harmonics) of the correlation function. 
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Uncertainty source 

*out (%) 

*ride (%) 

*lone (%) 

CF representation & 

Background parametrization 

5-32 

4-22 

4-35 

Fit-range dependence 

10 

8 

10 

7l + 7l + VS. 7T 7T 

3 

3 

3 

Momentum resolution correction 

3 

3 

3 

Two-track cut variation 

< 1 

< 1 

< 1 

Coulomb correction 

< 1 

< 1 

< 1 

Total correlated 

12-34 

9-24 

11-36 

Total 

12-34 

11-24 

12-36 


Table 2: List of contributions to the systematic uncertainty of the femtoscopic radii extracted via GGG fits. Values 
are averaged over kj and multiplicity except for the first row where a minimum-maximum range is shown. 


Uncertainty source 

*out(%) 

^de (%) 

tffon, (%) 

CF representation & 

Background parametrization 

4-18 

3-14 

8-20 

Fit-range dependence 

10 

6 

10 

71+71 + VS. 7T 7t 

3 

3 

3 

Momentum resolution correction 

3 

3 

3 

Two-track cut variation 

< 1 

< 1 

< 1 

Coulomb correction 

< 1 

< 1 

< 1 

Total correlated 

11-21 

7-16 

13-23 

Total 

12-21 

8-16 

14-23 


Table 3: List of contributions to the systematic uncertainty of the femtoscopic radii extracted via EGE fits. Values 
are averaged over kj and multiplicity except for the first row where a minimum-maximum range is shown. 


The same functional form for both of them was used for the fitting procedure. However, the implemen¬ 
tation of the fitting procedure is quite different: log-likelihood for Cartesian vs. regular % 2 fit for SH, 
3D Cartesian histogram vs. three ID histograms, cubic or spherical fitting range in the (<7out><7side,9iong) 
space. Therefore the fits to the two representations may react in a systematically different way to the 
variation of the fitting procedure (fit ranges, Bowler-Sinyukov approximation, etc.). 


The fitting procedure requires the knowledge of the non-femtoscopic background shape and ma gnit ude. 
Two models were used to estimate it, EPOS 3.076 liofill and PYTHIA 6.4 (Perugia-0 tune) 137, 38], as 
described in Sec. 13.21 


In addition the correlation function shape is not ideally described by a Gaussian form. The EGE form is 
better (lower % 2 values for the fit), but still not exactly accurate. As a result the fit values depend on the 
fitting range used in the procedure of radius extraction. We have performed fits with an upper limit of 
the fit range varied between 0.3 GeV/c and 1.1 GeV/c. 

The three effects mentioned above are the main sources of systematic uncertainty on the radii. Their 
influence, averaged over the event multiplicity and pair kj, is given in Tableland Tabled The back¬ 
ground parametrization and the CF representation effects lead to systematic uncertainties less than 10% 
at low kj and up to 35% for large kj and low multiplicities. In particular', the radius could not be reliably 
extracted for the two highest kj ranges in the lowest multiplicity range, therefore these two sets of radii 
are not shown. Moreover, radii obtained with the background parametrization from PYTHIA are always 
larger than the ones obtained with the EPOS parametrization. These uncertainties are correlated between 
kj ranges. Similarly, the radii from the narrow fit range are always on average 10% higher than the ones 


13 




















Two-pion femtoscopy in p-Pb collisions 


ALICE Collaboration 


from the wide fit range. This also gives a correlated contribution to the systematic uncertainty. The 
final radii are calculated as an average of four sets of radii - the two representations with both EPOS 
and PYTHIA background parametrization. The systematic uncertainties are symmetric and equal to the 
largest difference between the radius and one of the four sets of radii. 

The effect of the momentum resolution on the correlation function was studied using a Monte Carlo 
simulation. For tracks with a low pj, below 1 GeV/c, the momentum resolution in the TPC is better than 
1%. Smearing of the single particle momenta reduces the height and increases the width of the correlation 
function. It was estimated that this effect changes the reconstructed radius by 2% for a system size of 2 
fm and 3% for a size of 3 fm. Therefore, the 3% correlated contribution from momentum resolution is 
always added to the final systematic uncertainty estimation. 

Smaller sources of systematic uncertainties include those originating from the difference between posi¬ 
tively and negatively charged pion pairs (around 3%), track selection variation (less than 1%) and from 
the Coulomb factor (less than 1%). All the uncertainties are added in quadrature. 


4 Results 

4.1 Three-dimensional radii 


We have extracted R out , R s id e , and P| ong in intervals of multiplicity and kj, which results in 26 radii in 
each direction. The fit procedure did not allow to reliably extract values of the radii for the two highest 
kj ranges in the 60-90% multiplicity class. For the GGG fit, they are shown in Fig.|7] The radii are in the 
range of 0.6 to 2.4 fm in all directions and universally decrease with kj. The magnitude of this decrease 
is similar for all multiplicities in the out and long directions, and is visibly increasing with multiplicity 
in the side direction. The radii rise with event multiplicity. The plot also shows data from pp collisions 
at = 7 TeV ltl6ll at the highest multiplicities measured by ALICE, which is slightly higher than the 
multiplicity measured for the 20-40% V0A signal range in the p-Pb analysis. At small kj the pp radii are 
lower by 10% (for side) to 20% (for out ) than the p-Pb radii at the same (dY c h/dq) 1/3 . At high kj the 
difference in radius grows for /? out , while for /?i ong the radii for both systems become comparable. The 
distinct decrease of radii with kj is observed both in pp and p-Pb. 


The correlation strength X increases with kj from 0.44 to 0.58 for the collisions with highest multiplic¬ 
ities. It is also higher for low multiplicity collisions, with a difference of 0.1 between collisions with 
highest multiplicities and lowest multiplicities. A non-constant X parameter as a function of kj is an 
indication of a non-Gaussian shape of the correlation function. The correlation functions are normalized 
to the ratio of the number of pairs in the signal and background histograms. The positive correlation at 
low q has to be then compensated by the normalization parameter N, which is in the range of 0.9-1.0. 
The % 2 /ndf for the three-dimensional fit is on the order of 1.2. 


The extracted background parameters indicate that this contribution increases with kj and decreases with 
multiplicity, which is consistent with qualitative expectations for the mini-jet effect. The shape of the 
background is not spherical, leading to finite contributions to the (2,0) and (2,2) components. The 
constant shift in these components, given by (f) 1 and respectively, is only significant for the (2,0) 
component in lower multiplicities. 

The corresponding fit results for the EGE fit are shown in Fig. [8] In the side direction the radii are 
consistent with the GGG results. The radii in the out and long directions are not Gaussian widths in 
this case and cannot be directly compared to previous fits. However, all the trends are qualitatively the 
same in both cases: radii increase with event multiplicity and decrease with pair transverse momentum. 
The values are 10% (for side and long) to 20% (for out) higher than those measured in pp collisions at 
similar event multiplicity Q. The X parameter for the EGE fit is on the order of 0.7 for the SH fits, 
and growing from 0.7 at low kj to approximately 0.9 at the highest kj for the Cartesian fit, therefore 
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Fig. 7: (Color online) Femtoscopic radii (GGG fit) as a function of the pair transverse momentum kj for four mul¬ 
tiplicity classes. For comparison, radii from high-multiplicity pp collisions at yfs = 7 TeV 1150 and 4 predictions 
for p-Pb I24| 25] are shown as crosses and lines, respectively. Top, middle, and bottom panel show R ouU A M de. and 
R | ong radii, respectively. The points for multiplicity classes 20-40% and 40-60% have been slightly shifted in Ti¬ 
to r visibility. 


significantly higher than in the Gaussian case. The observation that A. is closer to unity when moving 
from GGG to EGE fits is expected, as the EGE fit describes the shape of the correlation much better at 
low q and therefore better accounts for the non-Gaussian tails in the source function. 


4.2 Model comparisons 


Hydrodynamic model calculations for p-Pb collisions 1124 25], shown as lines in Fig. [7] predict the ex¬ 
istence of a collectively expanding system. Both models employ two initial transverse size assumptions: 
A'mit = 1.5 fm and A' m ii = 0.9 fm, which correspond to two different scenarios of the energy deposition 
in the wounded nucleon model 12511 . The resulting charged particle multiplicity densities (dAch/dp) of 
45 [25] and 35 [240 are equal or higher than the one in the ALICE 0-20% multiplicity class. The calcula- 
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Fig. 8: (Color online) Femtoscopic radii (EGE fit) as a function of the pair transverse momentum k-\- for four 
multiplicity classes. For comparison, radii from high-multiplicity pp collisions at -y/s = 7 TeV lfl5l] are shown. 
Top, middle, and bottom panel show Rg Ul , Rf i( j e , Rf ong radii, respectively. The data points for the multiplicity 
classes 20-40% and 40-60% have been slightly shifted in kj for visibility. 


tions for R out overestimate the measured radii, while the ones with large initial size strongly overpredict 
the radii. The scenarios with lower initial size are closer to the data. For R s id e , the calculations are in 
good agreement with the data in the highest multiplicity class, both in magnitude and in the slope of the 
kj dependence. Only the Shapoval et al. [24] calculation for large initial size shows higher values than 
data. For f?i ong , calculations by Bozek and Broniowski [250 overshoot the measurement by at least 30% 
for the most central data, while those by Shapoval et al. are consistent within systematics. Again, the 
slope of the kj dependence is comparable. The study shows that the calculation with large initial size 
is disfavored by data. The calculations with lower initial size are closer to the experimental results, but 
are still overpredicting the overall magnitude of the radii by 10-30%. Further refinement of the initial 
conditions may lead to a better agreement of the models with the data, especially at large multiplicities. 
The slope of the kj dependence is usually interpreted as a signature of collectivity. Interestingly, it is 
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very similar in data and the models in all directions, which suggests that the system dynamics might be 
correctly modeled by hydrodynamics. 

Also in the data the source shape is distinctly non-Gaussian. Further studies would require examination 
of the source shape in p-Pb collision models to see if similar deviations from a Gaussian form are 
observed. 


The CGC approach has provided a qualitative statement on the initial size of the system in p-Pb col¬ 
lisions, suggesting that it is similar to that in pp collisions 11271. 281 1. The measured radii, at high mul¬ 
tiplicities and low kj, are 10-20% larger than those observed at similar multiplicities in pp data. For 
lower multiplicities the differences are smaller. These differences could still be accommodated in CGC 
calculations. Furthermore, the evolution of the slope of the kj dependence is similar between pp and 
p-Pb collisions in the side direction. Another similarity is the distinctly non-Gaussian shape of the 
source, which in pp and p-Pb is better described by an exponential-Gaussian-exponential form. It ap¬ 
pears that data in p-Pb collisions still exhibit strong similarities to results from pp collisions. However 
some deviations, which make the p-Pb more similar to A-A collisions, are also observed, especially 
at high multiplicity. The differences between small systems such as pp and p-Pb and peripheral A-A 
data are most naturally explained by the significantly different initial states in the two scenarios. Dedi¬ 
cated theoretical investigation of this issue is needed for a more definite answer, which may be able to 
accommodate both CGC and hydrodynamic picture. 


4.3 Comparison to the world data 


In Fig. [9]the results from this analysis of the p-Pb data from the LHC (red filled circles) are compared to 
the world heavy-ion data, including results obtained at lower collision energies, as well as to results from 
pp collisions from ALICE and STAR. It has been observed 111311 that the 3-dimensional femtoscopic radii 
scale roughly with the cube root of the measured charged-particle multiplicity density not only for a sin¬ 
gle energy and collision system, but also across many collision energies and initial system sizes. The pp 
and A-A datasets show significantly different scaling behavior, although both are linear in (dAch/dr)) 1 ^ 3 . 

The p-Pb radii agree with those in pp collisions at low multiplicities. With increasing multiplicity, the 
radii for the two systems start to diverge. An analysis of one-dimensional averaged radii in pp, p-Pb and 
Pb-Pb collisions using the 3-pion cumulant correlations technique reveals that the multiplicity scaling for 
p-Pb lies between pp and Pb-Pb trends [21], which is consistent with results presented here. On the other 
hand, the deviation of the correlation shape from Gaussian is similar to that observed in pp collisions, 
and unlike the shapes observed in A-A collisions. 


5 Conclusions 

We reported on the three-dimensional pion femtoscopic radii in p-Pb collisions at ^Ann = 5.02 TeV, 
measured in four multiplicity and seven pair momentum intervals. The radii are found to decrease with 
k'\ in all cases, similar to measurements in A-A and high-multiplicity pp collisions. The radii increase 
with event multiplicity. At low multiplicities they are comparable to the pp values, while at higher 
multiplicities and low pair transverse momentum they are larger by 10-20%. However, they do not reach 
the values observed in A-A collisions at lower energies. The high multiplicity data are compared to 
predictions from two models, both of them incorporating a fast hydrodynamic expansion of the created 
medium. They overpredict the values of the R out and /?i ong parameters, however the introduction of a 
smaller initial size results in a better description. The values of the P s ide parameter and the slope of 
the kj dependence of the radii are in reasonable agreement. The models based on the CGC formalism 
suggest sizes similar to those obtained in pp data. The observed differences of about 10-20% for high 
multiplicity p-Pb collisions might not exclude this scenario. The observed non-Gaussian shape of the 
correlation is also similar in the pp and p-Pb collision systems. 
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. (k T ) = 0.25 GeV/c 



0 5 10 

<d/v ch /dp) 1/3 



< dA/ ch /d / 7 > 1/3 



* STAR Au-Au s NN = 200 GeV 

* STAR Cu-Cu s^ N = 200 GeV 
STAR Au-Au s NN = 62 GeV 

0 STAR Cu-Cu \[s^ = 62 GeV 

* CERES Pb-Au \|^ JN = 17.2 GeV 

* ALICE Pb-Pb ^ = 2760 GeV 
■ ALICE pp is =7000 GeV 

□ ALICE pp fs = 900 GeV 
STAR pp is = 200 GeV 

* ALICE p-Pb \[s NN = 5020 GeV 


Fig. 9: (Color online) Comparison of femtoscopic radii (Gaussian), as a function of the measured charged-particle 
multiplicity density, measured for various collision systems and energies by CERES [42y, STAR [1451,1461 53], 
PHENIX J3, and ALICE [0]. 
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